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We explore the constant intensity cut method that is widely used for the derivation of the cosmic 
ray energy spectrum, for comparisons of data obtained at different atmospheric depths, for mea- 
suring average shower profiles, and for estimates of the proton-air cross section from extensive air 
shower data. The constant intensity cut method is based on the selection of air showers by charged 
particle or muon size and therefore is subject to intrinsic shower fluctuations. We demonstrate 
that, depending on the selection method, shower fluctuations can strongly influence the character- 
istics of the selected showers. Furthermore, a mixture of different primaries in the cosmic ray flux 
complicates the interpretation of measurements based on the method of constant intensity cuts. 
As an example we consider data published by the Akeno Collaboration. The interpretation of the 
Akeno measurements suggests that more than 60 — 70% of cosmic ray primaries in the energy range 
IQie - 10" eV are heavy nuclei. Our conclusions depend only weakly on the hadronic interaction 
model chosen to perform the simulations, namely SIBYLL and QGSjet. 
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I. INTRODUCTION 

Measuring extensive air showers (EAS) is currently the 
only way to study the cosmic ray spectrum and chemical 
composition at energies above 10^'* eV, as well as the 
basic properties of hadronic interactions at -^s above 1.8 
TeV. 

EAS can be detected with air shower arrays which 
measure densities of shower particles such as electrons, 
muons, photons, and sometimes hadrons arriving at the 
detector. These densities are typically fit to lateral dis- 
tribution functions to derive the total number of charged 
particles, electrons and muons A^^ at detector level. 
The particle numbers are functions of the primary cos- 
mic ray energy E and the mass number A of the primary 
particle, and depend on the atmospheric depth of the 
observation level. At energies E ^10^^ eV the shower 
evolution can also be directly observed by measuring the 
fluorescence light from the atmospheric nitrogen that is 
excited by the ionization of the charged shower particles. 
In the following we will concentrate on air shower arrays. 
Imaging methods such as fluorescence or Cherenkov light 
techniques will be discussed elsewhere [Q. 

One of the classical methods in the analysis of air 
shower data is the constant intensity cut method. The 
idea is based on the fact that, due to the isotropy of the 
primary cosmic ray flux, showers generated by primary 
particles of the same energy and composition will arrive 
at the detector with the same frequency, assuming 100% 
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detection efficiency. Selecting showers arriving at the 
detector with the same frequency under different zenith 
angles allows the measurement of the mean longitudinal 
shower profile. At large atmospheric depths after shower 
maximum, the shower size decreases approximately ex- 
ponentially with depth with a length scale commonly re- 
ferred to as the attenuation length. 

On the other hand, selecting showers with the same 
features (i.e. shower size, muon size etc.) at observation 
level and different incident angles allows the measure- 
ment of the absorption length which determines how the 
flux of the selected showers decreases with atmospheric 
depth. 

Measurements of the attenuation length are commonly 
used to correct observed particle densities to those of 
equivalent vertical showers. By unfolding the geometry- 
related attenuation of showers an experiment can use the 
measured intensities of showers with fixed size to derive 
the primary all-particle flux. 

The absorption length is inherently related to the mean 
free path of the EAS initiating primary particle. For 
example, the rate of proton air showers having the flrst 
interaction point (Xjnt) at a slant depth greater than 
X decreases as exp(— AT/Aint), where Aint is the mean 
free path for p — air collisions. On this basis, several 
methods of extracting the p — air cross section 23 from 
measurements of EAS have been applied in air shower 
experiments @ ||, ||, §. 

Air shower arrays cannot measure the depth of the 
flrst interaction of the primary particle generating the 
observed shower, which directly relates to the mean free 
path. The decrease with zenith angle of the frequency 
of showers having the same electron iVg and muon iV^ 
sizes at observation level is studied instead. In the ab- 
sence of intrinsic shower fluctuations these measurements 
would reflect the depth distribution of primary interac- 
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tions. However, the longitudinal development of showers 
is itself subject to large fluctuations. To disentangle these 
fluctuations from those of the flrst interaction point is not 
an easy task. 

This problem is usually addressed by introducing a co- 
efficient (fc) which relates the observed shower absorption 
length (Aobs) and the inelastic cross section through the 
equation Ai^t = kx Aobs @ • The numerical value of k has 
to be obtained from simulations of EAS. The coefficient k 
reflects the influence of the features of the hadronic inter- 
actions model on the fluctuations in the shower develop- 
ment. The value of k depends on the cross sections, sec- 
ondary particle multiplicity and elasticity in the hadronic 
interaction model. Due to the necessary extrapolation of 
hadronic multiparticle production to unmeasured regions 
of the phase space and to high energy, the extracted cross 
section becomes model dependent. 

Further difficulties in determining the inelastic p-air 
cross section from EAS measurements are related to ex- 
perimental uncertainties and limitations in the determi- 
nation of the development of air showers, and also to the 
fact that the cosmic ray flux might be "contaminated" 
with primaries heavier than protons which in principle 
tend to decrease the observed mean free path. 

In this article we shall study the importance of intrin- 
sic shower fluctuations for the experimentally observed 
attenuation and absorption lengths by considering two 
examples: 

(i) the reconstruction of the primary cosmic ray spectrum 
using charged particle shower sizes, and 

(ii) the measurement of the proton-air cross section, fol- 
lowing closely the method applied flrst by the Akeno 
group |§, ^ which we call the constant N^j^, method. 

In the process of (ii) we found that the Akeno observa- 
tions can best be understood if there is a large fraction of 
heavy nuclei present in the cosmic ray beam from 10-100 
PeV. 

The paper is structured as follows. In Section || we 
study the possible errors introduced in the derivation 
of the primary cosmic ray spectrum by shower fluctu- 
ations when the constant intensity method is applied. 
Section III consists of three parts. Part A summarizes 
the basics of the constant Ng — Nf^^ method. In part 
B we describe the predictions of this method for proton 
induced showers and in part C we discuss the more realis- 
tic situation of a mixed primary cosmic ray composition. 
Section IV concludes the paper. 



II. DERIVATION OF THE COSMIC RAY 
ENERGY SPECTRUM 

The constant intensity cut method has been used for 
studies of the primary cosmic ray spectrum for at least 
40 years. Let us give as an example the interpretation of 
the results of the BASJE air shower array at Mt. Cha- 
caltaya performed in 1965 Since Mt. Chacaltaya 
is at an altitude of 5,220 m above sea level (540 g/cm^ 



TABLE I: Results of the fit to the relation between shower 
energy and size at observation level for vertical showers at dif- 
ferent depths Xv and assuming different cosmic ray primaries. 
The last row corresponds to the energy-shower size relation 
obtained for a mixed composition assuming equal fractions 
of p. He, CNO and Fe. The function used to do the fit is 
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depth), the shower size distributions obtained with the 
constant intensity cut could be used to estimate the size 
of the showers at shower maximum. Under the assump- 
tion that the size at maximum is proportional to the 
primary energy E, this gave directly the primary energy 
spectrum within a constant which was estimated to be 
~ 2 GeV/particle at shower maximum. 

In this work we apply the constant intensity cut 
method in a different way, similar to what more con- 
temporary experiments do (see for instance pO[|). 

For illustrative purposes we flrst apply the method 
assuming the primary spectrum is composed of pure 
protons. We have simulated proton-induced showers at 
zenith angles 9 —0, 15, 30 and 45 degrees down to the 
altitude of the Akeno array, corresponding to a vertical 
depth of ATv = 920 g/cm^. Shower energies were drawn 
from an differential injection spectrum in the energy 
range between lO^*' and 10^^ eV. 

We used a hybrid air shower simulation program to 
generate large samples of showers in an efficient and fast 
manner ||ll|. The hybrid method consists of calculating 
shower observables by a direct simulation of the initial 
part of the shower, tracking all particles of energy above 
Ethi — 0.01 E. Presimulated showers for all subthreshold 
particles are then superimposed after their first interac- 
tion point is simulated. The subshowers are described 
with parametrizations that give the correct average be- 
havior, and at the same time describe the fluctuations in 
shower development of both electrons and muons. 

The procedure we use to reconstruct the primary spec- 
trum is the following. First, from the simulations we ob- 
tain the relation between shower energy (E) and shower 
size at observation level in vertical showers (N^)- This 
relation is given in Table |. Our simulations predict the 
shower size spectra at Akeno level for different zenith 
angles, and we treat them as if they were actual experi- 
mental data, replacing the detector induced fluctuations 
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FIG. 1: Energy resolution of proton-induced showers for dif- 
ferent zenith angles. The SIBYLL 2.1 hadronic generator code 
was used to simulate the showers. 

in logiQ Ne by a Gaussian resolution function of width 
AlogiQA^e — 0.05. Since we do not simulate showers 
of energy below 10^^ eV there is an artificial break in 
the size spectra at low energy. To avoid it we choose 
Ne = 10^ as a threshold value above which our "array" 
is fully efficient, and we only deal with showers having 
Nf, above this value. 

We apply cuts at constant shower intensity, and by 
studying the decrease with zenith angle of the size cor- 
responding to each intensity, we obtain the shower at- 
tenuation length, which we use to estimate the shower 
size at = 0° from the known size at zenith angle 9. 
(We checked that the attenuation length obtained from 
the simulated data in this way agrees with the attenua- 
tion length of the averaged profile of the input showers.) 
This is the classical integral application of the method. 
Given we can use the previously obtained relation 
between E and to estimate the energy of each indi- 
vidual shower. The energy spectrum for different zenith 
angles can then be reconstructed and compared to the 
injected spectrum. 

Fig. |l| shows the shower energy resolution we achieve 
with this procedure. The distribution of differences be- 
tween reconstructed and injected energies is highly asym- 
metric. The asymmetry depends on zenith angle as show- 
ers at larger zenith angles are further away from their 
maximum where the fluctuations in shower size are small- 
est. There is a clear tendency to misreconstruct showers 
of a certain injected energy assigning them a higher en- 
ergy. 

In Fig. y we compare the reconstructed and injected 
spectra for different zenith angles. The spectra are mul- 
tiplied by E"^-^ for a better resolution. Although we draw 
showers from an E~^ differential spectrum, the cut in Nf, 
decreases the contribution of the lower energy cosmic rays 
and creates the turnaround seen in Fig. |^. The larger 
the zenith angle, the higher the shower energy must be 
to exceed the N^, cut, producing a strongly zenith angle 



FIG. 2: Reconstructed (points) and injected (histograms) en- 
ergy spectra (multiplied by E'^'^ at different zenith angles. 
Note that larger zenith angles contribute to the derived cos- 
mic ray spectrum in a limited energy range because of the 
Ne > 10^ threshold. 

dependent energy cut. 

Besides, as mentioned previously, with increasing 
zenith angle the shower is sampled further away from 
shower maximum and the fluctuations in Ne grow, pro- 
ducing a broadening of the injected spectrum. These two 
effects, however, vanish almost completely in the recon- 
structed spectra. 

Imposing a cut in iVe produces a corresponding cut 
in E through the relation between energy and size at 
observation level. In our case Ne > 10^ implies that 
logio E > f6.36, 16.66, 17.16 at zenith angle 0, 30 and 45 
deg respectively. As a consequence all the low energy in- 
jected events are reconstructed with energies above these 
values, and they pile-up rendering a simple power law. A 
fit to the reconstructed spectra reveals that the differen- 
tial spectral index decreases by only about 2 — 3%. This 
is an important effect which should be present in exper- 
iments reconstructing the spectrum which use a relation 
between shower energy and shower size at observation 
depth, for instance at a certain distance from the core of 
the shower. 

Although the spectral shape is preserved almost com- 
pletely there is a slight, but noticeable difference in the 
absolute normalization of the spectra reconstructed from 
showers at different zenith angles. The energies derived 
from showers at non vertical angle are always overesti- 
mated, which leeds in principle to an artificially increased 
normalization. Fig. |^ shows that increase for an angle of 
30°, however the normalization from the 45° showers is 
again lower. This may be due to the energy resolution 
distribution that peaks well below for these showers. 
Independently of the exact reason for the changing nor- 
malization. Fig. ^ demonstrates that the exact energy 
derivation depends strongly on the shape of the shower 
fluctuations. Since these fluctuations change with the 
atmospheric depth, so does the reconstruction accuracy. 
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As a whole, though, the method works quite well. 

The reconstruction of the shower energy from Nf. is af- 
fected even less by the shower fluctuations if the depth 
of the detector is close to the depth of shower maximum. 
In this case the danger is in the inclusion of showers 
that have not yet achieved their maximum development. 
Such showers may introduce a significant bias when their 
size is converted to vertical size by using the attenuation 
length. In the example discussed above, the mean depth 
of shower maximum is ~ 650 g/cm^ with a standard devi- 
ation of ~ 70 g/cm^, so that only a small fraction ^ 0.5% 
of the vertical showers have their maxima below obser- 
vation level. This fraction is much smaller for inclined 
showers. 

Further difficulties in the reconstruction arise in the 
more realistic case when the primary spectrum consists 
of a mixed composition of different nuclei. The heavier 
the primary nucleus, the further away is the observation 
level from shower maximum. From this point of view the 
spectrum reconstruction for a mixed primary composi- 
tion is analogous to using proton showers to very large 
zenith angles. 

To explore this point we have simulated a primary com- 
position consisting of equal fractions of protons. He, CNO 
and Fe. We obtained the E — Nf. relation at 920 g/cm^ 
(shown in table ^ from this particular mixture of nuclei, 
and applied the same procedure as before to reconstruct 
the primary spectrum. The reconstruction is again af- 
fected by the cut in Nf. as explained above. The result 
is that both the spectral index and the normalization of 
the reconstructed spectrum differ from the corresponding 
values in the injection spectrum by only a few percent. 

It is important to note that we have made use of our 
prior knowledge of the injected primary composition to 
reconstruct the spectrum. We expect the reconstruction 
method to work either when the primary composition is 
known or when a composition independent energy esti- 
mator is used. (The density at 600 m from the cores of 
large showers is an example of a measure of shower en- 
ergy chosen because of its relative insensitivity to primary 
mass Q.) We have checked that using the E — Ne re- 
lation obtained for pure protons tends to underestimate 
the normalization of a mixed spectrum. The reason is 
that a shower of energy E initiated by a heavy nucleus 
has on average a smaller size at observation level than a 
proton shower of the same energy. A smaller energy is 
then assigned to the shower when the E — relation 
for protons is used. As expected, we observed the op- 
posite behavior when the E — Nf, relation for pure iron 
(also shown in table ||) is used to reconstruct the mixed 
composition spectrum. 

The most difficult case is obviously that of changing 
chemical composition. Because of the changes of the 
spectrum normalization for different primary nuclei the 
shape of the spectrum can also be derived incorrectly. 
The composition and spectrum then have to be recon- 
structed simultaneously from different shower parame- 
ters. 



We conclude that the integral application of the con- 
stant intensity cuts method for the derivation of the pri- 
mary cosmic ray spectrum works well when the cosmic 
ray composition is known. The use of wrong composition 
models can lead to erroneous conclusions for the energy 
spectrum, mostly in the determination of its normaliza- 
tion. The method is not strongly affected by the intrinsic 
fluctuations in the shower in the energy range we have 
explored. 



III. 



THE CONSTANT A^,, - A^^ METHOD 



The total column density of atmosphere available for 
shower development increases with the incident angle 9 
as secO. The total number of electrons at a fixed slant 
depth after the shower maximum reflects the stage of 
evolution of the shower. If shower fluctuations were ab- 
sent, selecting showers of fixed energy at different zenith 
angles which have the same electron size A^e j would a pri- 
ori guarantee that they have developed through the same 
column of atmosphere between the first interaction point 
and observation level. The selected showers would only 
differ in the depth at which the first primary p-air inter- 
action had occurred. The proton-air interaction length 
Xp-air and the corresponding cross section Up-air would 
then be measured. The fact is, however, that Nf. does 
have large fluctuations. 

To address this problem and select showers of fixed pri- 
mary energy, experiments often require that the showers 
have the same muon size at observation level N^. Unlike 
electrons, the number of muons A^^ remains almost con- 
stant after maximum, and hence it is a good estimator of 
the primary energy at essentially any observation depth 
below shower maximum. Selecting the showers with large 
electron sizes within the same A'^ bin increases the prob- 
ability that they are induced by protons. The higher 
the size is, the lower is the contamination from heavier 
primaries. 

Once showers are selected in this way, the frequency 
(/) of showers falling in a given (A"^, A'e) bin is measured 
for different zenith angles. The ratio of the frequency 
of selected showers at two zenith angles (6*1 and 62) is 
related to the observed absorption length by 



f{N^,N,,e^) 



2) 



Av 



exp — - — — (sec 6*1 — sec 6*2 ) 

(1) 

where is the vertical depth of the detector. 

Clearly, intrinsic fluctuations of the shower profile will 
change the relation between the first interaction point 
and the electron and muon number at larger depths. In 
the following we will study how such fluctuations influ- 
ence the observed absorption length, using detailed, up- 
to-date hadronic interactions models. For definiteness we 
will concentrate on the implementation of the method as 
used by the Akeno group 1 1, ^ . A similar procedure was 
also used more recently by the EAS-TOP Collaboration 
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A. Application to proton-induced showers 

In general, the primary cosmic ray flux consists of nu- 
clei of a variety of mass numbers. Since we want to study 
the constant — N^^ method itself we simplify the prob- 
lem and start with the assumption that all primary par- 
ticles are protons. If the procedure does not give the 
correct cross section for a purely proton flux, the correct 
derivation for a mixed cosmic ray composition would be 
impossible. 

We have performed simulations of proton-induced 
showers at several zenith angles and calculated the fre- 
quency of showers having iV^ muons and iVe electrons at 
observation level. The detector was chosen to be located 
at Akeno altitude, corresponding to a vertical depth of 
Xy — 920 g/cm^. Shower energies were drawn from an 
E^^ differential spectrum in the energy range between 
10^^ and 10^^ eV. We performed the simulations for fixed 
zenith angles {9 =0, 15, 30 and 45 deg.). We thus sim- 
plify the problem once again by neglecting the errors in- 
troduced by the experimental shower zenith angle recon- 
struction. 

To study the dependence on the hadronic interaction 
model we have performed our simulations with two mod- 
els, namely SIBYLL 2.1 |l|, |l§ and QGSjet98 |l|. The 
two models give similar predictions for the shower devel- 
opment in the energy range 10^^ — 10^^ eV [0. As will 
be shown, our results depend only weakly on the choice 
of model. 

For a realistic simulation of the observed shower pa- 
rameters one should account for the experimental uncer- 
tainty and fluctuations due to the detector. A detailed 
simulation of the biases and efficiencies of the detectors 
is beyond the scope of this paper. It requires the use of 
specifically designed Monte Carlo programs of each par- 
ticular ground array. We replace instead the detector 
induced fluctuations in log^Q and log^g by Gaus- 
sian resolution functions of widths Alogj^g -^m = 0-1 a-nd 
Alogj^QA^e = 0.05 respectively. These are the experi- 
mental errors reported by the Akeno group Q| . For each 
of the simulated showers modified log^Q and log^o 
are sampled according to the theoretical values and the 
detector resolution. 

It is possible and even likely that, due to different en- 
ergy thresholds and absorbing materials, the experimen- 
tal definition of TVe does not coincide exactly with the 
corresponding quantity that the Monte Carlo generates. 
In such a case one should apply a correction factor to 
achieve a full reproduction of the experimental result. 
As we show further below, however, a possible small dis- 
crepancy in the definition of iVg would not alter the con- 
clusions of this study. 

We have simulated a sample of 500,000 proton show- 
ers, comparable to the statistics of the full event sample 
reported by the Akeno collaboration in B . 
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FIG. 3: Top panel: Simulated energy distribution of proton- 
initiated showers after applying the muon cut. The selected 
showers have logjp A^^ between 5.25 and 5.45. The muons 
have energy above Ef"^ = 1 GeV x sec 6 at 920 g/cm^ The 
distribution is shown for different zenith angles. The SIBYLL 
2.1 hadronic generator code was used to simulate the showers. 
Bottom panel: Same as top panel after applying the constant 
intensity cut. 



We apply the constant A'^e — method by first se- 
lecting showers which have log^o between 5.25 and 
5.45 at observation level, the muon number of the first 
bin of the Akeno analysis [||. Only muons with energy 
E^> 1 GeV X sec 61, which is the muon energy threshold 
of the Akeno experiment Q], are considered. In the top 
panel of Fig. ^ we plot the energy distribution of showers 
selected after applying the muon cut. For fixed primary 
energy the mean muon number is smaller at large zenith 
angles as compared to vertical showers. This is due both 
to the dependence of the energy threshold on zenith an- 
gle, and to the increase of the probability for muon de- 
cay when the zenith angle increases and muons have to 
traverse more column depth of atmosphere. As a conse- 
quence, the energy distribution shows a dependence on 
zenith angle, i.e. selecting showers with the same number 
of muons does not perfectly guarantee that they have the 
same energy distribution. One can correct for the iV^ at- 
tenuation by the constant intensity cut method shifting 
slightly the log^g -^/^ bin so that the intensity of showers 
is the same at all zenith angles. This is equivalent to a 
correction of the shower muon longitudinal profile as a 
function of the zenith angle. The bottom panel of Fig. ^ 
shows the energy distribution of the selected showers af- 
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ter applying the constant intensity cut method. The con- 
stant intensity cut method works almost perfectly, giving 
a distribution of selected shower energies independent of 
zenith angle. The width of the shower energy distribu- 
tion is determined by the width of the bin, and by 
the shower to shower fluctuations. 

Once showers of the same energy have been selected, 
we select in addition showers with constant Ng, as was 
done by the Akeno team. The top panel of Fig. |^ shows 
the A^e spectra of showers having log^p -^m between 5.25 
and 5.45 for the four nominal zenith angles. The two ver- 
tical lines mark the bin in logj^Q A'e chosen by the Akeno 
collaboration to perform their analysis of the showers in 
the 5.25-5.45 bin in log^g N^^. The bottom panel of Fig. | 
shows the energy distribution of showers in the selected 
Ne bin. It shows that the energy estimate that was very 
good after the Nf^ bin selection with constant intensity 
cuts, is now again angle dependent. 

We have compared the Nf,{E) dependence calculated 
by the Monte Carlo code to the experimental one used 
by the Akeno experiment ||l^. In the log^Q A^e bin 6.8 
- 7.0 the Akeno formula gives log]^g(£'/eV) = 16.40 and 
the Monte Carlo code using SIBYLL 2.1 yields 16.43. 
Both values are estimated at the center of the bin. 
The differences when using the N^{E) are higher. The 
Akeno collaboration [|l6| derived \ogi^^{E / eV) — 16.25 
in the logig A^^, bin 5.25 - 5.45. SIBYLL 2.1 gives 



logio(i;/eV) = 16.55 and QGSjet98 logio(i;/eV) = 16.48 
in the same log^Q A'^ bin. 

In Fig. ^ we further illustrate this dependence by show- 
ing the frequency ratios of the showers in Fig. ^, in de- 
pendence of the selected electron size. The ratio is only 
plotted for adjacent zenith angles. It depends strongly on 
the A'e bin used for shower selection. According to Eq. (|l|) 
this ratio should be constant over a certain range in A^e 
for all different zenith angle combinations. Fig. || shows 
that the A^e ranges where the ratio is approximately con- 
stant depend on the shower angle combination. For both 

The bin in cho- 
sen by Akeno for the cross section analysis is clearly in a 
region where the intensity ratios depend strongly on A"e. 

The ratio of two ratios R for different combinations of 
zenith angles can be used as a consistency check of the 
results. The expected values of the double ratio are (see 
Eq. (0)) 



SIBYLL 2.1 and QGSjet98 models they lie above log^o 
of 7.4 for the logio bin 5.25-5.45 



logioi?(15°,30°) _ sec(30) - sec(15) 



logioi?(0°,15° 



sec(15) — sec(O) 



logio i?(3G°, 45°) sec(45) - sec(30) 



3.4 



7.4 



(2) 



(3) 



logio R{0°, 15°) sec(15) - sec(O) 

It is easy to verify that this is not the case for the se- 
lected showers falling in the bin in log^g A'e between 6.8 
and 7.0. These numbers are, however, the approximate 
scaling factors when comparing the plateaus of the ratios 
in Fig. ||. The figure suggests that zenith angle depen- 
dent bins in electron size should be used in order to get 




FIG. 4: Top panel: Shower size at 920 g/cm^ depth in 
proton-initiated showers having between 10^'^^ and 10^'*^ 
muons at 920 g/cm^ The size distribution is shown for show- 
ers initiated at different zenith angles. Histograms correspond 
to showers simulated using SIBYLL 2.1, and points to showers 
simulated with QGSjet98. Bottom panel: Energy distribution 
of the showers falling in the A^e bin indicated by the vertical 
bars in the top panel. 



an angular-independent value of Aobs , provided of course 
the selected showers are initiated by protons as in our 
simulation. 

A similar consistency check can be performed by plot- 
ting the observed shower intensity as a function of sec 6. 
Fig. ^ shows the intensity of proton-initiated showers 
falling in the \ogiQ Nf^ = 5.25 — 5.45 bin for different 
A^e bins. The deviation from straight lines, which are 
expected for exponential attenuation of showers with 
sec 6*, demonstrates that the constant A'e — A'^ method 
fails to select similar showers, unless a large value of 
A'e is selected. The intensity of the selected show- 
ers certainly does not decrease as exp(— X/Aobs) in the 
(logio logio ^e) = (5.25 - 5.45, 6.8 - 7.0) bin. Formal 
fits of the three higher A'e bins plotted in Fig. ||give Aobs 
values of 138 ±36, 85 ±8 and 68 ±3 g/cm^ for log^g A^e of 
7.0-7.2, 7.2-7.4 and 7.4-7.6 respectively. Compared to the 
proton- air cross section of 456 mb in SIBYLL 2.1 these 
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FIG. 5: Ratios of number of proton-initiated showers having 
between 10^'^^ and lO''"^ muons and electron size A^e at 920 
g/cm^ as a function of N^- Histograms correspond to showers 
simulated using SIBYLL 2.1, and points to showers simulated 
with QGSjet98. 
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FIG. 6: Zenith angle dependence of the intensity of proton- 
induced showers having constant logj^Q A^^ = 5.25 — 5.45 and 
constant logjQ A^e for different values of log^gNe- Empty 
squares logj^Q A^e = 6.8—7.0, filled squares logjQ A'^e = 7.0 — 7.2, 
empty circles logj^Q A'^^ = 7.2 — 7.4 and filled circles logjQ A^e ~ 
7.4 - 7.6. Showers were simulated with SIBYLL 2.1. The 
points are joined by straight lines to guide the eye. To avoid 
overlapping, the results for different bins were multiplied 
by different arbitrary factors. 



values lead to k values of 2.60 ± 0.67, 1.61 ± 0.14, and 
1.28 ±0.06. 

The analogous analysis carried out with QGSjet98 for 
the same muon number bin and logj^Q = 7.4 — 7.6 
gives Aobs = 69 ± 2 which corresponds to a fc-factor of 
1.26 ±0.04. Within the statistical uncertainty this value 
agrees with the one derived from SIBYLL 2.1 simula- 



FIG. 7: Distribution in Xobs — Xint of the showers that fall in 
the (logio Ar^,logio Afe)=(5.25-5.45,6.8-7.0) bin. 



tions. The weak model dependence is not unexpected. 
The energy of the showers considered here is only one 
order of magnitude higher than the equivalent energy of 
the Tevatron collider. Both models were tuned to re- 
produce the Tevatron measurements and predict rather 
similar muon and electron numbers for E ^ 10^^ eV. 



B. Shower fluctuations 

The ultimate reason why the constant — method 
does not work is that the discussed shower selection is 
dominated by the intrinsic fluctuations in shower de- 
velopment. This is illustrated in Fig. ^ in which we 
plot the distribution of "shower lengths" of showers with 
(logio A^M'logio^e) = (5-25 - 5.45,6.8 - 7.0) for differ- 
ent angles. We arbitrarily define the shower length as 
the difference between the slant depth of observation 
level (Xobs) and the slant depth of the first interaction 
point. If the longitudinal shower profile were not biased 
by the selection criteria all four histograms would be very 
similar. They are instead very different and demostrate 
that the selection is on the width of shower development 
rather than on the depth of the first interaction X-mt- 
Indeed, for all angles most of the selected showers have 
their first interaction point near the top of the atmo- 
sphere. Even in the larger Ne bins the situation is not 
qualitatively different, as can be seen in Fig. ||. 

In table |] we give the average values and the widths of 
the distributions shown in Fig 0. For the method not to 
be dominated by intrinsic fluctuations, the average value 
of Xobs— ^int should be independent of zenith-angle. The 
increase of atmospheric depth from 6* = 0° to 45° should 
lead to a shift of the first interaction point by about 400 
g/cm^. Due to the fluctuations the actual mean shift is 
only by about 110 g/cm^. 

A way to quantify which part of the longitudinal 
shower development, namely the first few interactions 
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or the latest interactions, contributes most to the fluc- 
tuation in shower length, is calculating the average val- 



10000 



ues and widths of the distributions of 



Xs„t and 



Xobs — ^max- This is shown in columns 3 and 4 in the 
table. The tail of the shower contributes more to the 
overall fluctuation in shower length than the first few 
interactions, although the contribution depends on the 
zenith angle. In terms of the ratio of cr/AX the angular 
dependence is bigger for (Xobs — -''^max), where it changes 
from 0.054 for vertical showers to 0.136 for showers de- 
veloping under 45°. 



C. Composition 

In contrast to our findings summarized in Fig. |6[ the 
Akeno Collaboration reports a secO dependence of the 
observed frequencies of showers selected by the constant 
A^e ~ method which is compatible with an exponen- 
tial attenuation (see Fig. 1 in Ref. Q). By looking at 
Figs. ^ and ^ it is clear that the intensity of proton show- 
ers in the nominal bin is not large enough at small zenith 
angles to produce a straight line. Proton showers pen- 
etrate too much in the atmosphere and thus have large 
electron sizes at observation level. In principle, this state- 
ment depends on the hadronic interaction model used in 
the simulations. A model which predicts the same muon 
number at lower shower energy can lead to an increase of 
the vertical shower intensity in the considered bin. How- 
ever, simulations with QGSjet, which predicts the largest 
muon multiplicity among the contemporary hadronic in- 
teraction models, show that this conclusion is unchanged 
if the muon multiplicity in the considered energy range 
is increased by up to 20%. 

A way to increase the intensity of the selected show- 
ers would be to "contaminate" the sample with heavy 
primaries. These give rise to showers which are less 
penetrating, shifting the distribution in electron size to 
smaller values. This is illustrated in the top panel of 
Fig. ^ for a primary composition consisting of 85% Fe, 
10% CNO, 4% He and 1% protons. The bottom panel 
of Fig. ^ shows the contributions to the total electron 
size distribution from showers initiated by the different 
primaries. It is remarkable, and to our understanding a 
coincidence, that the tail of the total distribution in elec- 
tron size has roughly the same slope as the tail of the 
contribution from proton-induced showers alone. 

Fig. ^ shows the ratio of the electron size distributions 
shown in the top panel of Fig. ^ for adjacent zenith an- 
gles. In the Nf. range where the ratios are fiat, they have 
numerical values very similar to the expected ratios for 
protons shown in Fig. ^ in the corresponding plateau re- 
gions, i.e. although the primary spectrum is dominated 
by heavy primaries the analysis method gives a value of 
the cross section similar to, but somewhat higher than 
that obtained for pure protons. 

The composition we chose in this analysis is com- 
pletely "ad-hoc" . In particular we have assumed an 




FIG. 8: Top panel; Same as Fig. ^ for a primary cosmic 
ray composition consisting of 85% Fe, 10% CNO, 4% He and 
1% protons. The bottom panel shows the 6 — 30° size dis- 
tribution that is plotted in the top panel illustrating how the 
different cosmic ray primaries contribute to build it up. 




FIG. 9: Ratios of number of proton-initiated showers having 
between 10^'^^ and lO'^ *'^ muons and electron size A'^e at 920 
g/cm^ as a function of Ne for a primary cosmic ray spectrum 
consisting on 85% Fe, 10% CNO, 4% He and 1% protons. 
Showers were simulated with SIBYLL 2.1. 



energy-independent ratio of the different elemental con- 
tributions. However, we have repeated the analysis for 
many different combinations of primary fractions, and 
only those with a large fraction of iron produce an expo- 
nential decrease of the intensity of showers with zenith 
angle. This is shown in Fig. O. 
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TABLE II: Average values and standard deviations (in parenthesis) of the distributions of Xobs — A^int (see also Fig. 
Amax — Aint and Aobs — Amax for proton-initiated showers belonging to the (logj^Q Apjlogj^Q Ae) = (5.25-5.45,6.8-7.0) bin. 



Q [deg.] 


Aobs [g/cm^] 


(Aobs - Aint) [g/cm^] 


(Amax - Aint) [g/cHl^] 


(Aobs - Amax) [g/cm^] 





920.0 


881.3 (35.7) 


581.4 (31.1) 


300.1 (42.6) 


15 


952.2 


911.1 (37.6) 


585.8 (31.6) 


325.2 (42.3) 


30 


1062.3 


1002.0 (50.6) 


613.7 (39.5) 


388.2 (51.3) 


45 


1301.1 


1152.9 (109.7) 


699.7 (95.1) 


453.8 (92.5) 



100 



O.I 



0.01 



1% p, 4% He, 10% CNO, 85% Fe 
25% p, 25% He, 25% CNO, 25% Fe 
85% p, 10% He, 4% CNO, 1% Fe 
100% p 



IV. 



DISCUSSION AND CONCLUSIONS 
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FIG. 10: Zenith angle dependence of the intensity of show- 
ers having constant (logj^Q A^^, logj^Q iVe)=(5.25-5.45,6.8-7.0) 
for different primary compositions. Showers were simulated 
with SIBYLL 2.1. The points are joined by straight lines to 
guide the eye. To avoid overlapping of the results for differ- 
ent compositions, they are multiplied by different arbitrary 
factors. 



We have not attempted to perform a fit to the intensity 
versus zenith angle using the different fractions of pri- 
maries as parameters in the fit. Such an analysis would 
require the use of a true detector Monte Carlo simulation 
and could only be performed by the experimental group. 
However from the combinations we have experimented 
with we conclude that at least 60-70% of iron is needed 
to produce a straight line in the nominal N iVgbin. The 
Akeno measurements then strongly indicate that a large 
fraction of iron is present in the cosmic ray spectrum in 
the energy region between IQi^-lO" eV. This resuh IS m 
qualitative agreement with recent analyses of the region 
around and above the knee in the cosmic ray spectrum 
(for example, KASCADE and HiRes-MIA |l| mea- 
surements, see also [^). Finally, we note that an early 
analysis using the method of constant intensity cuts 
reached a similar conclusion about heavy composition in 
this energy range based on data from the BASJE air 
shower experiment on Mt. Chacaltaya |pt|. 



We have investigated the influence of air shower fluctu- 
ations on the widely used constant intensity cut method. 
We consider two types of applications: the classic inte- 
gral approach to the derivation of the cosmic ray energy 
spectrum, and the differential 7V^,iVe cut used for the 
derivation of the proton-air cross section. 

We find that the constant intensity cuts method can 
work for comparisons of data taken at different atmo- 
spheric depths and different angles. This is however pos- 
sible only when the chemical composition of the primary 
cosmic rays is well known. The use of incorrect chemical 
composition can lead to a shift in the normalization of the 
energy spectrum. In the case of energy dependent com- 
position the normalization errors for different cosmic ray 
fiux components could also affect significantly the derived 
spectral index 7. Such shifts are also possible close to the 
detector iVg threshold, where measurements at different 
zenith angles would detect showers of different composi- 
tion. The larger the zenith angle 0, the lighter would be 
the composition of the detected showers. 

The infiuence of shower fluctuations is much bigger 
when the constant intensity cut is used in a differential 
way to compare showers with same electron and muon 
sizes detected at different zenith angles. The selection by 
the muon size N p_ with constant intensity cuts is indeed a 
very good method and leads to a good angle independent 
energy selection. This is the result of the much slower 
absorbtion of the shower muons as well as to the smaller 
iVp fluctuations in showers with fixed primary energy. 

The constant Nf.—N^ method, which is used for deriva- 
tion of the proton-air production cross section, is domi- 
nated by fluctuations even in the case of a pure proton 
composition. The accuracy of this method improves with 
the selection of showers with large Nf. for a fixed bin, 
where an experiment would run out of statistics. A pos- 
sible improvement of the method would be to use Monte 
Carlo shower simulations to determine zenith angle de- 
pendent bins. This, however, would represent a new 
method which is very different from the original idea of 
constant intensity cuts. 

The Akeno data that were used for the derivation of the 
proton-air production cross section can be interpreted in 
terms of cosmic ray composition. The angle independent 
exponential slope of the shower absorption length indi- 
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cates a substantial fraction of heavy primaries in the en- 
ergy range of 10^^ - 10^^ eV. The comparison to showers 
simulated with the QGSjet98 and SIBYLL 2.1 hadronic 
interactions models require 60-70% of iron in the primary 
cosmic ray flux to explain the absorption length derived 
by the Akeno group. 

These conclusions depend only mildly on the hadronic 
interaction model used in the simulation. 
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